\documentclass{ctexart}
\usepackage[affil-it]{authblk}
\usepackage[dvipsnames]{xcolor}
\usepackage{amssymb}
\usepackage{verbatim}
\usepackage{indentfirst}
\usepackage{hyperref}
\usepackage{amsmath}
\usepackage{geometry}
\usepackage{xeCJK}
\usepackage{float}
\usepackage{listings}
\usepackage{graphicx} % Required for figures
\usepackage{subcaption} % Required for subfigures
\geometry{margin=1.5cm, vmargin={0pt,1cm}}
\setlength{\topmargin}{-1cm}
\setlength{\paperheight}{29.7cm}
\setlength{\textheight}{25.3cm}

\begin{document}
% ==================================================
\title{Numerical Analysis Homework \#3}

\author{周川迪 Zhou Chuandi 3220101409}
\affil{强基数学2201}

\date{\today}

\maketitle

\begin{abstract}
  3.4.1 Theoretical questions

  Complete question 3.4.1 I - VIII on page 31 of the lecture notes in English. 

  Theorems or Corollaries are referred from \textit{handoutsNUMPDEs-2024-08-17}
\end{abstract}



\section*{I. Determine whether a natural spline}
For \[ \ s(x) = \begin{cases} 
p(x) & \text{if } x \in [0, 1], \\
q(x) = (2 - x)^3 & \text{if } x \in [1, 2].
\end{cases} \]

There is\[q'(x)=-3(x-2)^2,\ q''(x)=6(2-x)\]

To satisfy \(s \in \mathbb{S}_2^3 \), We have 
\[ p(0)=0,\ p(1)=q(1)=1,\ p'(1)=q'(1)=-3,\ p''(1)=q''(1)=6\]

Use Hermite interpolation to determine \(p\)
\[
\begin{array}{c|cccc}
  0 & 0 & & & \\
  1 & 1 & 1 & & \\
  1 & 1 & -3 & -4 & \\
  1 & 1 & -3 & 3 & 7 
\end{array}
\]

We get \[p(x)=0+1(x-0)-4(x-0)(x-1)+7(x-0)(x-1)(x-1)=7x^3-18x^2+12x\]
\[\Rightarrow s''(0)=p''(0)=-36 \neq 0,\ s''(1)=q''(1)=6 \neq 0\]

Thus, \(s\) is not a natural spline.



\section*{II. Quadratic spline interpolation}
\subsection*{(a)}
Set \(p_i=s|_{[x_i,x_{i+1}]},\ i=1,\cdots,n-1 \), each \(p_i\) need \(3\) conditions to be uniquely determined.

But currently we only have \(3n-4\) conditions: \[p_i(x_i)=f_i,\ p_i(x_{i+1})=f_{i+1},\ i=1,\cdots,n-1\] \[p'_i(x_{i+1})=p'_{i+1}(x_{i+1}),\ i=1,\cdots,n-2\]

Thus, we need an additional condition to determine \(s\) uniquely.


\subsection*{(b)}

\[
\begin{array}{c|ccc}
  x_i & f_i & & \\
  x_i & f_i & m_i & \\
  x_{i+1} & f_{i+1} & \frac{f_{i+1}-f_i}{x_{i+1}-x_i} & \frac{f_{i+1}-f_i-m_i(x_{i+1}-x_i)}{(x_{i+1}-x_i)^2} 
\end{array}
\]

We get\[p_i(x) = f_i + m_i(x-x_i)+\frac{f_{i+1}-f_i-m_i(x_{i+1}-x_i)}{(x_{i+1}-x_i)^2}(x-x_i)^2\]

\subsection*{(c)}
When \(m_i\) is given, \(p_i\) can be uniquely determined. Then we can obtain \(m_{i+1}\):
\[m_{i+1}=s'(x_{i+1})=p_i(x_{i+1}) = 2\frac{f_{i+1}-f_i}{x_{i+1}-x_i}-m_i\]

Since \(m_1\) is given, we can compute \(m_2,m_3,\cdots,m_{n-1}\) iteratively.



\section*{III. Constructing a natural cubic spline}
For \( s_1(x) = 1 + c(x + 1)^3 \), there is \[ s'_1(x) = 3c(x + 1)^2,\ s''_1(x) = 6c(x + 1) \]

Set \(s_2=Ax^3+Bx^2+Cx+D\), then \[s'_2(x)=3Ax^2+2Bx+C,\ s''_2(x)=6Ax+2B\]

\(s\) is a natural cubic spline and knot \(-1\) has been OK by \(s_1\).

For the other two knots, \(s_2\) should satisfy that:

\[
\begin{cases} 
  s_2(0)=s_1(0)=1+c &\Rightarrow D=1+c \\
  s'_2(0)=s'_1(0)=3c &\Rightarrow C=3c \\
  s''_2(0)=s''_1(0)=6c &\Rightarrow B=3c \\
  s_2(1)=s(1)=-1 &\Rightarrow A+B+C+D=-1 \\
  s''_2(1)=0 &\Rightarrow 3A+B=0 
\end{cases}
\]

We solve out that \(c=-\frac{1}{3}\).



\section*{IV. Interpolating cosine function with natural cubic spline}
\subsection*{(a)}
For \( f(x) = \cos\left(\frac{\pi}{2}x\right) \), there is \(f(-1)=0,\ f(0)=1,\ f(1)=0\).

Set the natural cubic spline as \[s(x) = \begin{cases}
  s_1(x)=a_1x^3+b_1x^2+c_1x+d_1, & x \in [-1, 0], \\
  s_2(x)=a_2x^3+b_2x^2+c_2x+d_2, & x \in [0, 1]. 
\end{cases} \]

There are eight variables we have eight following equations to give the unique spline:
\begin{align*}
  &s_1(-1)=f(-1)=0,\ s_1(0)=f(0)=1,\ s''_1(-1)=0 \\
  &s_2(0)=f(0)=1,\ s_2(1)=f(1)=0,\ s''_2(1)=0 \\
  &s'_1(0)=s'_2(0),\ s''_1(0)=s''_2(0) 
\end{align*}

Solving it by hand is too tiring, so we can observe that \(s_1'(0)\) and \(s_2'(0)\) might be \(0\) and guess
\[s(x) = \begin{cases}
  s_1(x)=a_1x^3+b_1x^2+1, & x \in [-1, 0], \\
  s_2(x)=a_2x^3+b_2x^2+1, & x \in [0, 1]. 
\end{cases} \]

We soon find that here is the right solution:
\[s(x) = \begin{cases}
  s_1(x)=-\frac{1}{2}x^3-\frac{3}{2}x^2+1, & x \in [-1, 0], \\
  s_2(x)=\frac{1}{2}x^3-\frac{3}{2}x^2+1, & x \in [0, 1]. 
\end{cases} \]

This gives the natural cubic spline interpolant to f on those knots.


\subsection*{(b)}
We have \( a=-1,\ b=1\). And for \(s\), \[s''(x) = \begin{cases}
  -3x-3, & x \in [-1, 0], \\
  3x-3, & x \in [0, 1]. 
\end{cases} \]
\[\Rightarrow \int_a^b [s''(x)]^2 \, dx = \int_{-1}^0 (-3x-3)^2 \, dx + \int_0^1 (3x-3)^2 \, dx = 6\]

\subsubsection*{(i)}
Use Newton's interpolation to obtain the quadratic polynomial \(g\):
\[
\begin{array}{c|ccc}
  -1 & 0 & &  \\
  0 & 1 & 1 & \\
  1 & 0 & -1 & -1
\end{array}
\]

We get \[g(x)=x+1-1(x+1)x=-x^2+1,\ g''(x)=-2\]

\[\Rightarrow \int_{-1}^{1} [g''(x)]^2 \, dx = 8 > \int_a^b [s''(x)]^2 \, dx\]

\subsubsection*{(ii)}
\[f(x) = \cos\left(\frac{\pi}{2}x\right),\ f''(x) = -\frac{\pi^2}{4}\cos\left(\frac{\pi}{2}x\right)\]

Then,
\[\int_{-1}^{1} [f''(x)]^2 \, dx = \frac{\pi^4}{16}\int_{-1}^{1} \cos^2\left(\frac{\pi}{2}x\right)\, dx = \frac{\pi^4}{16}\int_{-1}^{1} \frac{1+\cos(\pi x)}{2}\, dx = \frac{\pi^4}{16} = \approx 6.09 \]

Thus,
\[\int_{-1}^{1} [f''(x)]^2\, dx > \int_a^b [s''(x)]^2 \, dx\]



\section*{V. Analysis of quadratic B-spline}

\subsection*{(a)}
By \textbf{Definition 3.23},
\[ B_i^0(x) = \begin{cases} 
  1 & \text{if } x \in (t_{i-1}, t_i], \\
  0 & \text{otherwise}.
\end{cases} \]

\[B_i^1(x) = \frac{x - t_{i-1}}{t_{i} - t_{i-1}} B_i^0(x) + \frac{t_{i+1} - x}{t_{i+1} - t_i} B_{i+1}^0(x) = \begin{cases} 
  \frac{x-t_{i-1}}{t_i-t_{i-1}} & x \in (t_{i-1}, t_i], \\
  \frac{t_{i+1}-x}{t_{i+1}-t_i} & x \in (t_i, t_{i+1}], \\
  0 & \text{otherwise}. 
\end{cases} \]

By \textbf{Definition 3.21},
\[B_i^1(x) = \hat{B}_i(x) \]

Examining the cases where \(x\) lies within the intervals respectively, the result is:
\[B_i^2(x) = \frac{x - t_{i-1}}{t_{i+1} - t_{i-1}} B_i^1(x) + \frac{t_{i+2} - x}{t_{i+2} - t_i} B_{i+1}^1(x) = \begin{cases} 
\frac{(x-t_{i-1})^2}{(t_{i+1}-t_{i-1})(t_i-t_{i-1})}, & x \in (t_{i-1}, t_i], \\
\frac{(x-t_{i-1})(t_{i+1}-x)}{(t_{i+1}-t_{i-1})(t_{i+1}-t_i)} + \frac{(t_{i+2}-x)(x-t_i)}{(t_{i+2}-t_i)(t_{i+1}-t_i)}, & x \in (t_i, t_{i+1}], \\
\frac{(t_{i+2}-x)^2}{(t_{i+2}-t_i)(t_{i+2}-t_{i+1})}, & x \in (t_{i+1}, t_{i+2}], \\
0, & \text{otherwise}.
\end{cases} \]
 

\subsection*{(b)}
Set \(b_1,\ b_2,\ b_3\) as follows:
\[B_i^2(x) = \begin{cases} 
b_1(x) = \frac{(x-t_{i-1})^2}{(t_{i+1}-t_{i-1})(t_i-t_{i-1})}, & x \in (t_{i-1}, t_i], \\
b_2(x) = \frac{(x-t_{i-1})(t_{i+1}-x)}{(t_{i+1}-t_{i-1})(t_{i+1}-t_i)} + \frac{(t_{i+2}-x)(x-t_i)}{(t_{i+2}-t_i)(t_{i+1}-t_i)}, & x \in (t_i, t_{i+1}], \\
b_3(x) = \frac{(t_{i+2}-x)^2}{(t_{i+2}-t_i)(t_{i+2}-t_{i+1})}, & x \in (t_{i+1}, t_{i+2}], \\
0, & \text{otherwise}.
\end{cases} \]

We Calculate \(b'_2(x)\) first: \begin{align*}
  b'_2(x) &= 
  -\frac{x-t_{i-1}}{(t_{i+1}-t_{i-1})(t_{i+1}-t_{i})}
  +\frac{t_{i+1}-x}{(t_{i+1}-t_{i-1})(t_{i+1}-t_{i})}
  -\frac{x-t_{i}}{(t_{i+2}-t_{i})(t_{i+1}-t_{i})}
  +\frac{t_{i+2}-x}{(t_{i+2}-t_{i})(t_{i+1}-t_{i})} \\
  &=\frac{t_{i+1}+t_{i-1}-2x}{(t_{i+1}-t_{i-1})(t_{i+1}-t_{i})}
  + \frac{t_{i+2}+t_{i}-2x}{(t_{i+2}-t_{i})(t_{i+1}-t_{i})}
\end{align*}


Calculate the derivatives: \begin{align*}
  & b'_1(t_i) = \frac{2}{t_{i+1}-t_{i-1}} \\
  & b'_2(t_i) = \frac{t_{i+1}+t_{i-1}-2t_i}{(t_{i+1}-t_{i-1})(t_{i+1}-t_{i})} + \frac{1}{t_{i+1}-t_{i}} = \frac{2}{t_{i+1}-t_{i-1}} \\
  & b'_2(t_{i+1}) = -\frac{1}{t_{i+1}-t_{i}} + \frac{t_{i+2}+t_{i}-2t_{i+1}}{(t_{i+2}-t_{i})(t_{i+1}-t_{i})} = -\frac{2}{t_{i+2}-t_{i}} \\
  & b'_3(t_{i+1}) = -\frac{2}{t_{i+2}-t_{i}}
\end{align*}

Since \(b'_1(t_i) = b'_2(t_i),\ b'_2(t_{i+1})=b'_3(t_{i+1}) \), \( \frac{\mathrm{d}}{\mathrm{d}x}B_i^2(x) \) is continuous at \(t_i\) and \(t_{i+1} \)


\subsection*{(c)}

\(\forall x \in (t_{i-1},t_{i}],\ \frac{\mathrm{d}}{\mathrm{d}x}B_i^2(x) = b'_1(x) = \frac{2(x-t_{i-1})}{(t_{i+1}-t_{i-1})(t_{i}-t_{i-1})} > 0 \)

\(\forall x \in (t_{i},t_{i+1}],\ \frac{\mathrm{d}}{\mathrm{d}x}B_i^2(x) = b'_2(x) \)

\( b'_2(x) \) is linear, with \( b'_2(t_i) > 0 \) and \( b'_2(t_{i+1}) < 0 \) 

Thus, there exists unique \(x^* \in (t_{i-1},t_{i+1})\ s.t.\ \frac{\mathrm{d}}{\mathrm{d}x}B_i^2(x^*) = 0 \)

By the expression of \(b'_2(x)\), we reached that
\[x^*= \frac{t_{i+2}t_{i+1}-t_{i}t_{i-1}}{t_{i+2}+t_{i+1}-t_{i}-t_{i-1}} \]

\subsection*{(d)}
\(\forall x \in [t_{i+1},t_{i+2}),\ \frac{\mathrm{d}}{\mathrm{d}x}B_i^2(x) = b'_3(x) = \frac{-2(t_{i+2}-x)}{(t_{i+2}-t_{i})(t_{i+2}-t_{i+1})} < 0 \)

Since \(B^2_i(x)=0,\ \text{otherwise}.\)

Thus, fix any \(i\), \(\forall x \in \mathbb{R}\), we have:
\[ B^2_i(x) \geq 0 \]
\[ B^2_i(x) \leq B^2_i(x^*) = b_2(x^*)\]

It can be calculated that \[b_2(x^*) = \frac{t_{i+2} - t_{i-1}}{t_{i+2} - t_{i-1} + t_{i+1} - t_i } < 1 \]

Thus, \(B^2_i(x) \in [0,1) \).

\subsection*{(e)}
Here is the python code for plotting: (Not too much so I put it here.)
\begin{lstlisting}[   % 进行参数设置
  language=Python, % 设置语言
  basicstyle=\ttfamily, % 设置字体族
  breaklines=true, % 自动换行
  keywordstyle=\bfseries\color{NavyBlue}, % 设置关键字为粗体，颜色为 NavyBlue
  morekeywords={}, % 设置更多的关键字，用逗号分隔
  emph={self}, % 指定强调词，如果有多个，用逗号隔开
     emphstyle=\bfseries\color{Rhodamine}, % 强调词样式设置
     commentstyle=\itshape\color{black!50!white}, % 设置注释样式，斜体，浅灰色
     stringstyle=\bfseries\color{PineGreen!90!black}, % 设置字符串样式
     columns=flexible,
     numbers=left, % 显示行号在左边
     numbersep=2em, % 设置行号的具体位置
     numberstyle=\footnotesize, % 缩小行号
     frame=single, % 边框
     framesep=1em % 设置边框与代码的距离
 ]  
  import numpy as np  
  import matplotlib.pyplot as plt  

  def Bspline(n, t, i, x): # spline function
    if (n == 0):
      return (x > t[i - 1]) * (x <= i) * 1.0
    else:
      return (x - t[i - 1]) / n * Bspline(n - 1, t, i, x) + (t[i] + n - x) / n * Bspline(n - 1, t ,i + 1, x)
  n = 2 # degree of splines
  m = 6 # number of splines (range of i)
  t = np.zeros(100) # knots
  for i in range(0, 90): # set knots
    t[i]=i
  x = np.linspace(0, m+3, 1000)
  Z = np.linspace(0, m+3, m+4) 
  plt.plot(Z, np.zeros(len(Z)), 'o')
  for i in range(2, 1+m):
    plt.plot(x, Bspline(n, t, i, x), label=f'$B_{i}^{n}$')
  plt.legend()
  plt.savefig("V.png")
\end{lstlisting}

\begin{figure}[H]
  \centering
  \includegraphics[width=0.8\textwidth]{"V.png"}
  \caption{The graph of \(B_i^2(x)\) for \(t_i = i,\ i = 2, 3, \dots, 6\)}
\end{figure}



\section*{VI. Verification of Theorem 3.32 for n=2}
What we aim to verify is:
\[ (t_{i+2} - t_{i-1})[t_{i-1}, t_i, t_{i+1}, t_{i+2}](t - x)^2 _+ = B_i^2(x) \]

Here \begin{align*}
  \text{LHS}&=[t_{i},t_{i+1},t_{i+2}](t-x)^2 _+ - [t_{i-1},t_{i},t_{i+1}](t-x)^2 _+ \\
  &=\frac{\frac{(t_{i+2}-x)^2 _+ - (t_{i+1}-x)^2 _+}{t_{i+2}-t_{i+1}} - \frac{(t_{i+1}-x)^2 _+ - (t_{i}-x)^2 _+}{t_{i+1}-t_{i}} }{t_{i+2}-t_{i}} 
   -\frac{\frac{(t_{i+1}-x)^2 _+ - (t_{i}-x)^2 _+}{t_{i+1}-t_{i}} - \frac{(t_{i}-x)^2 _+ - (t_{i-1}-x)^2 _+}{t_{i}-t_{i-1}} }{t_{i+1}-t_{i-1}}
\end{align*}

It can be calculated that
\[
\begin{cases}
  x \in (-\infty,t_{i-1}],\ & \text{LHS}=0 \\
  x \in (t_{i-1},t_{i}],\ & \text{LHS}=\frac{(x-t_{i-1})^2}{(t_{i+1}-t_{i-1})(t_i-t_{i-1})}\\
  x \in (t_{i},t_{i+1}],\ & \text{LHS}=\frac{(x-t_{i-1})(t_{i+1}-x)}{(t_{i+1}-t_{i-1})(t_{i+1}-t_i)} + \frac{(t_{i+2}-x)(x-t_i)}{(t_{i+2}-t_i)(t_{i+1}-t_i)}\\
  x \in (t_{i+1},t_{i+2}],\ & \text{LHS}=\frac{(t_{i+2}-x)^2}{(t_{i+2}-t_i)(t_{i+2}-t_{i+1})}\\
  x \in (t_{i+2},+\infty),\ & \text{LHS}=0\\
\end{cases}
\]

By the expression of \(B_i^2(x) \) we have reached in \textbf{Problem V}, LHS=RHS. Done.




\section*{VII. Scaled integral of B-splines}
The support of \(B_i^0(x)\) is \( (t_{i-1},t_i] \).

By \textbf{Theorem 3.23}, \[B_i^{n+1}(x) = \frac{x - t_{i-1}}{t_{i+n} - t_{i-1}} B_i^n(x) + \frac{t_{i+n+1} - x}{t_{i+n+1} - t_i} B_{i+1}^n(x)\]

By induction, the support of \(B_i^n(x)\) is \( (t_{i-1},t_{i+n}] \). 

So the scaled integral of \( B_i^n(x) \) is
\[I_i^n= \int_{t_{i-1}}^{t_{i+n}} \frac{B_i^n(x)}{t_{i+n} - t_{i-1}} \, dx \]

By \textbf{Theorem 3.34} (Derivatives of B-splines), take \(n+1\):
\[\frac{\mathrm{d}}{\mathrm{d}x} B_i^{n+1}(x) = \frac{(n+1) B_i^{n}(x)}{t_{i+n} - t_{i-1}} - \frac{(n+1) B_{i+1}^{n}(x)}{t_{i+n+1} - t_i}\]

Calculate the difference between scaled integrals:
\begin{align*}
  I_{i+1}^n-I_i^n &= \int_{t_{i-1}}^{t_{i+n}} \frac{B_i^n(x)}{t_{i+n} - t_{i-1}} \, dx -\int_{t_{i}}^{t_{i+n+1}} \frac{B_{i+1}^n(x)}{t_{i+n+1} - t_{i}} \, dx \\
  &= \int_{t_{i-1}}^{t_{i+n+1}} \frac{B_i^n(x)}{t_{i+n} - t_{i-1}} - \frac{B_{i+1}^n(x)}{t_{i+n+1} - t_{i}} \, dx \\
  &= \frac{1}{n+1} \int_{t_{i-1}}^{t_{i+n+1}} \frac{\mathrm{d}}{\mathrm{d}x} B_i^{n+1}(x) \, dx \\
  &=\frac{B_i^{n+1}(t_{i+n+1}) - B_i^{n+1}(t_{i-1})}{n+1} \\
  &=0
\end{align*}

\[\Rightarrow I_{i+1}^n=I_i^n,\ \forall i. \]

Therefore, the scaled integral of \( B_i^n(x) \) is independent of its index \(i\).



\section*{VIII. Symmetric Polynomials as Divided Differences}
The theorem we aim to prove is \textbf{Theorem 3.46}:
\begin{equation}\label{1}
  \tau_{m-n}(x_i, \cdots, x_{i+n}) = [x_i, \cdots, x_{i+n}] x^m
\end{equation}

By \textbf{Theorem 3.38},
\[ \tau_{m-n}(x_{i}, \cdots, x_{i+n}) = \sum_{i \leq i_1 \leq \cdots \leq i_{m-n} \leq i+n} x_{i_1} x_{i_2} \cdots x_{i_{m-n}}.\]

\subsection*{(a)}
Calculate the table of divided differences for \(m=4\) and \(n=2\).
\[
\begin{array}{c|ccc}
    x_i & x_i^4 & & \\
    x_{i+1} & x_{i+1}^4 & x_{i+1}^3 + x_{i+1}^2x_i + x_{i+1}x_i^2 + x_i^3 & \\
    x_{i+2} & x_{i+2}^4 & x_{i+2}^3 + x_{i+2}^2x_{i+1} + x_{i+2}x_{i+1}^2 + x_{i+1}^3 & x_{i}^2+ x_{i+1}^2+x_{i+2}^2 + x_ix_{i+1}+ x_ix_{i+2} + x_{i+1}x_{i+2}
\end{array}
\]

Thus, \[\tau_{2}(x_{i},x_{i+1},x_{i+2}) = x_{i}^2+ x_{i+1}^2+x_{i+2}^2 + x_ix_{i+1}+ x_ix_{i+2} + x_{i+1}x_{i+2} = [x_i,x_{i+1},x_{i+2}]x^4\]

Verified.

\subsection*{(b)}
Fix \(m \in \mathbb{N}^+\), assume \eqref{1} holds for some \(0 \leq n<m\), then by the lemma on the recursive relation,
\begin{align*}
  &(x_{i+n+1}-x_i)\tau_{m-n-1}(x_i, \cdots, x_{i+n+1}) \\
  =\ &
  \tau_{m-n}(x_i, \cdots, x_{i+n+1})-\tau_{m-n}(x_i, \cdots, x_{i+n})-
  \tau_{m-n}(x_i, \cdots, x_{i+n+1})+\tau_{m-n}(x_{i+1}, \cdots, x_{i+n+1})\\
  =\ & \tau_{m-n}(x_{i+1}, \cdots, x_{i+n+1})-\tau_{m-n}(x_i, \cdots, x_{i+n}) \\
  =\ & [x_{i+1}, \cdots, x_{i+n+1}] x^m - [x_i, \cdots, x_{i+n}] x^m
\end{align*}

We obtain that
\[\tau_{m-n-1}(x_i, \cdots, x_{i+n+1})=[x_i, \cdots, x_{i+n+1}]x^m\]

It implies that \eqref{1} also holds for \(n+1\).

Since \(\forall m \in \mathbb{N}^+\), when \(n=0\), there is \(\tau_m(x_i)=x_i^m=[x_i]x^m\), which satisfies \eqref{1}.

Therefore, by induction, \eqref{1} holds for all \( m \in \mathbb{N}^+\) and \(n=0,1,\cdots,m\).

Q.E.D.
\end{document}